#ifndef __PASE_SOL_H__
#define __PASE_SOL_H__

#include "pase_mg.h"
#include "pase_param.h"
#include "ops_eig_sol_gcg.h"
#include "ops.h"
#include "petscmat.h"
#include "pase_ops.h"

typedef struct PASE_MG_SOLVER_PRIVATE_
{
    //======================= multigrid =====================
    MULTIGRID_TYPE multigrid_type;
    PASE_MULTIGRID multigrid;
    int num_levels;
    // for GMG
    BOUNDARY_TYPE boundary_type;
    bool boundary_delete;
    int boundary_num;
    int *boundary_index;
    //======================= problem setting =====================
    int user_nev;       // 用户要算多少个特征值
    int pase_nev;
    int num_given_eigs; // 给了多少个特征值（目前都是0）
    //======================= PASE =====================
    PASE_TYPE solver_type;
    PC_TYPE pc_type;
    bool pc_eigendcp;
    int aux_coarse_level;
    int aux_fine_level;
    int max_cycle_count;
    //======================= 分组 =====================
    bool if_groups;
    int group_num;
    int group_id;
    MPI_Group group;
    MPI_Comm group_comm;
    int real_nev;           // 本组所计算的特征对的个数
    int nev_start, nev_end; // 本组所计算的特征对的起始编号
    int actual_start, actual_end;
    //======================= 分批 =====================
    bool if_batches;
    int batch_num;
    int batch_size;
    int more_batch_size;
    int real_batch_size;
    double rtol;
    int more_aux_nev;
    // change during solving
    int current_batch;
    int current_nev_start, current_nev_end;
    int current_cycle;
    int conv_nev;
    int nlock_auxmat_A;
    int nlock_auxmat_B;
    bool *convergence_check;
    //======================= smoothing =====================
    SMOOTHING_TYPE smoothing_type;
    int max_pre_count, max_post_count;
    int nlock_smooth;
    //======================= initial GCGE =====================
    int initial_level;
    int max_initial_direct_count;
    double initial_atol;
    double initial_rtol;
    //======================= aux GCGE =====================
    int max_direct_count;
    double current_shift;
    double aux_atol;
    double aux_rtol;
    double max_aux_rtol;
    double aux_rtol_coef;
    double *aux_eigenvalues;
    // double *init_eigenvalues;
    void *shift_factorization;
    //======================= eigen pairs =======================
    double *eigenvalues;
    void ***solution;
    //======================= tmp space =====================
    void ***cg_rhs;
    void ***cg_res;
    void ***cg_p;
    void ***cg_w;
    double *double_tmp;
    int *int_tmp;
    double *double_tmp_hrr;
    // for aux GCGE
    PASE_MultiVector aux_gcg_mv[4];
    // PASE_MultiVector rr_mv;
    // for precondition
    double *pre_double;
    void **pre_bh;
    //======================= ops ===========================
    OPS *gcge_ops;
    PASE_OPS *pase_ops;
    OPS *pase_ops_to_gcge;
    //======================= aux problem =====================
    PASE_MultiVector aux_sol;
    PASE_Matrix aux_A;
    PASE_Matrix aux_B;
    //======================= timer =====================
    double solver_setup_time;
    double get_initvec_time;
    double smooth_time;
    double build_aux_time;
    double prolong_time;
    double error_estimate_time;
    double aux_direct_solve_time;
    double precondition_setup_time;
    double precondition_solve_time;
    double total_solve_time;
    double total_time;
    //======================================================
    int print_level;
    bool if_error_estimate;
} PASE_MG_SOLVER_PRIVATE;
typedef PASE_MG_SOLVER_PRIVATE *PASE_MG_SOLVER;

int PASE_EigenSolver(double **eval, void ***evec, PASE_PARAMETER param);
int PASE_EigenSolve(PASE_MG_SOLVER solver , PASE_PARAMETER param, PASE_OPS *pase_ops);
int PASE_Mg_set_up(PASE_MG_SOLVER solver, PASE_PARAMETER param);
int PASE_Mg_solve(PASE_MG_SOLVER solver);
int PASE_Direct_solve(PASE_MG_SOLVER solver);
int PASE_Mg_cycle(PASE_MG_SOLVER solver);
int PASE_Mg_cycle_simple(PASE_MG_SOLVER solver);
int PASE_Mg_smoothing(PASE_MG_SOLVER solver, int times);
int PASE_Aux_direct_solve(PASE_MG_SOLVER solver);
int PASE_Aux_direct_solve_simple(PASE_MG_SOLVER solver);
int PASE_Mg_error_estimate(PASE_MG_SOLVER solver);

int PASE_Mg_find_shift(PASE_MG_SOLVER solver, int rr);
int PASE_Mg_prolong_from_Solution(PASE_MG_SOLVER solver, int object_level);
int PASE_Mg_choose_eigenpairs(PASE_MG_SOLVER solver);
int PASE_Mg_choose_eigenpairs_simple(PASE_MG_SOLVER solver);
int PASE_Mg_prolong_from_pase_aux_Solution(PASE_MG_SOLVER solver, int object_level);
int PASE_Mg_set_pase_aux_vector(PASE_MG_SOLVER solver);
int PASE_Mg_set_pase_aux_matrix(PASE_MG_SOLVER solver);
int PASE_Mg_set_pase_aux_matrix_simple(PASE_MG_SOLVER solver);

PASE_MG_SOLVER PASE_Mg_solver_create(PASE_PARAMETER param, PASE_OPS *pase_ops);
int PASE_Mg_pase_aux_vector_create(PASE_MG_SOLVER solver, int aux_vecnum, int vecnum);
int PASE_Mg_pase_aux_matrix_create(PASE_MG_SOLVER solver, int aux_vecnum);
int PASE_Mg_pase_aux_matrix_factorization(PASE_MG_SOLVER solver, PASE_PARAMETER param);

int PASE_Matrix_destroy_sub(PASE_Matrix *aux_A);
int PASE_Matrix_destroy_sub_wofactorization(PASE_Matrix *aux_A);
int PASE_MultiVector_destroy_sub(PASE_MultiVector *aux_sol);
int PASE_Mg_solver_destroy(PASE_MG_SOLVER solver);
int PASE_Mg_solver_destroy_wofactorization(PASE_MG_SOLVER solver, bool A, bool B);

int PASE_Mg_print_param(PASE_MG_SOLVER solver);
int PASE_Mg_print_timer(PASE_MG_SOLVER solver);

int PASE_DIRECT_GCGE(void *A, void *B, int flag, int nev, double atol, double rtol, int argc, char *argv[]);
int PASE_DIRECT_EPS(void *A, void *B, int flag, int nev, double tol, int argc, char *argv[]);

int PASE_SETUP(PASE_MG_SOLVER* solver, PASE_PARAMETER param);
int PASE_EIGENSOLVE(PASE_MG_SOLVER solver);
int PASE_DESTROY(PASE_PARAMETER *param,PASE_MG_SOLVER solver);
#endif
